library(Seurat)
library(Matrix)
library(tidyverse)
library(viridis)

Load the Seurat object

Before we load the object, we’re going to switck the rownames of our data to be the symbols. There is a 1-1 mapping of gene ID to symbol in these samples so we can do this with no loss of information.

dat = readRDS(file.path("data", "counts.rds"))
md = readRDS(file.path("data", "metadata.rds"))
g2s = readr::read_csv(file.path("metadata", "tx2gs.csv"),
                      col_names=c("gene", "symbol", "transcript")) %>%
  dplyr::select(-transcript) %>%
  unique()
rownames(dat) = g2s$symbol[match(rownames(dat), g2s$gene)]
seurat = CreateSeuratObject(dat, meta.data=md)

Overview

First we’ll look at the uncorrected data; we’ll just correct for mitochondrial percentage and nothing else. Then we’ll run PCA, clustering and look at the UMAP plots to get an idea of how much run-to-run and technology variation we are dealing with.

seurat = SCTransform(seurat, vars.to.regress = "pctMito", verbose = FALSE)
seurat = RunPCA(seurat, verbose = FALSE)
seurat = RunUMAP(seurat, dims = 1:30, verbose = FALSE)
p1 = DimPlot(seurat, group.by = "technology")
p1

A lot! We can see that the 10x and the inDrop data are mostly separate from each other, but there is some overlaps. We can see looking at the four samples that the four samples are very different from each other as well.

DimPlot(seurat, group.by = "sample")

So we are going to have some work to do trying to overlay these on top of each other. There are a few ways to do that. One is to just regress out the effect of the technology, the other is to do something more complex like CCA. We do need to do something about it though if we want to combine these two datasets. However, even without doing anything about the technological variation we can take a peek at the data and see if we can pick out some patterns.

uncorrected clustering

Here we’ll do a rough clustering and then look at the expression of our markers, without correcting for which technology or experiment they are from, to orient outselves.

clustering

ElbowPlot(seurat, ndims=50)

seurat = FindNeighbors(seurat, dims=1:25)
seurat = FindClusters(seurat, resolution=0.5)
## Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck
## 
## Number of nodes: 10605
## Number of edges: 332518
## 
## Running Louvain algorithm...
## Maximum modularity in 10 random starts: 0.9401
## Number of communities: 25
## Elapsed time: 0 seconds
DimPlot(seurat)

cell type markers

markers = readr::read_csv("metadata/cell_type_markers.csv")
knitr::kable(markers)
cell_identity genes
ISC Dl
ISC Smvt
ISC CG10006
pISC sna
ISCR2R5 polo
ISC/EB esg
EB E(spl)m3-HLH
EB E(spl)malpha-BFM
EB E(spl)mbeta-HLH
EB klu
EE pros
EE 7B2
EE_subtype AstA
EE_subtype AstC
EE_subtype CCHa1
EE_subtype CCHa2
EE_subtype Tk
EE_subtype NPF
EE_subtype Orcokinin
EE_subtype Dh31
aEC nub
aEC Myo31DF
mEC nub
mEC Myo31DF
pEC nub
pEC Myo31DF
aEC betaTry
pEC lambdaTry
mEC lab
mEC Vha100-4
cardia Pgant4
copper PGRP-LA
copper PGRP-LC
copper Apoltp
VM Vn
LFC Ilp3
LFC PGRP-SC1a
LFC PGRP-SC1b
LFC PGRP-SD
iron PGRP-SC2
iron ZIP1(Zip42C.1)

Below we can see that cluster 3 and cluster 21 are likely ISC+EB or EB cells. Cluster 9, 17 and 20 are liekly some subtypes of EE cells. Groups 1, 4, 11, 14 and 23 might be aEC cells. Groups 5, 13 and 15 might be pEC cells. Group 8 looks like mECs. And group 15 looks like cardia cells.

Cluster 12 looks like it might be LFC cells. Cluster 7 might be iron or copper cells. Cluster 6 looks like it could be an EC cell, but we don’t know what compartment.

We’re just missing 2, 10, 18, 19, and 22. Not bad.

cluster.averages = AverageExpression(seurat, return.seurat=TRUE, add.ident="sample")
DoHeatmap(cluster.averages, features=head(markers, 20)$genes, size=2)

DoHeatmap(cluster.averages, features=tail(markers, 21)$genes, size=2)

clusters = seurat@meta.data$seurat_clusters
rough_celltype = case_when(
  clusters %in% c(3, 21) ~ "ISC/EB",
  clusters %in% c(9, 17, 20) ~ "EE",
  clusters %in% c(0, 1, 4, 11, 14, 23) ~ "aEC",
  clusters %in% c(5, 13, 16) ~ "pEC",
  clusters %in% c(8, 24) ~ "mEC",
  clusters %in% c(6) ~ "EC",
  clusters %in% c(15) ~ "cardia",
  clusters %in% c(12) ~ "LFC",
  clusters %in% c(7) ~ "iron/copper",
  TRUE ~ as.character(clusters))
seurat@meta.data$rough_celltype = rough_celltype
Idents(seurat) = seurat@meta.data$rough_celltype
roughcelltype.averages = AverageExpression(seurat, return.seurat=TRUE, add.ident="sample")
DoHeatmap(roughcelltype.averages, features=head(markers, 20)$genes, size=2)

DoHeatmap(roughcelltype.averages, features=tail(markers, 21)$genes, size=2)

p3 = DimPlot(seurat, group.by="rough_celltype", label=TRUE)
p3

It looks like we are able to identify the cells across both technologies, which is very good.

knitr::kable(seurat@meta.data %>%
             group_by(rough_celltype, technology) %>%
             summarise(n=n()))
rough_celltype technology n
10 10x 76
10 inDrop 323
18 10x 31
18 inDrop 214
19 inDrop 243
2 10x 100
2 inDrop 683
22 10x 10
22 inDrop 77
aEC 10x 913
aEC inDrop 3142
cardia 10x 32
cardia inDrop 225
EC 10x 186
EC inDrop 316
EE 10x 306
EE inDrop 551
iron/copper 10x 206
iron/copper inDrop 266
ISC/EB 10x 141
ISC/EB inDrop 643
LFC 10x 102
LFC inDrop 221
mEC 10x 252
mEC inDrop 226
pEC 10x 624
pEC inDrop 496

unmarked clusters

This leaves us with clusters 2, 10, 18, 19 and 22 as unknown clusters.

seurat@meta.data %>%
  filter(is.na(rough_celltype)) %>%
  pull(seurat_clusters) %>%
  unique()
## factor(0)
## 25 Levels: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 ... 24

trypsin genes

tryp = c("deltaTry", "gammaTry", "alphaTry", "betaTry",
         "epsilonTry", "thetaTry", "kappaTry", "lambdaTry",
         "iotaTry", "etaTry", "zetaTry")
r1 = c("deltaTry", "gammaTry", "Muc68D")
r2 = c("deltaTry", "gammaTry", "alphaTry", "betaTry", "epsilonTry", "thetaTry")
r3 = c("thetaTry")
r4 = c("kappaTry", "lambdaTry", "iotaTry", "etaTry", "zetaTry")
r5 = c("iotaTry", "etaTry", "zetaTry")
crop = c("Cyp312a1", "Cyp4e3", "Spn27A", "spz", "pot")
compartments = unique(c(r1, r2, r3, r4, r5, crop))
FeaturePlot(seurat, features=r2)

FeaturePlot(seurat, features=r4)

FeaturePlot(seurat, features=r5)

FeaturePlot(seurat, features=crop)

DoHeatmap(cluster.averages, features=compartments, size=2)

It looks like group 12 and 15 might be from R1. Group 0, 1, 4, 11, 14, 19, 23 and 24 are from either R2 or R3. Group 7 looks like it might be from R3. Groups 5, 6 and 16 look like they are from R4. Group 18 looks like it is from the crop.

clusters = seurat@meta.data$seurat_clusters
compartment = case_when(
  clusters %in% c(12, 15) ~ "R1",
  clusters %in% c(0, 1, 4, 11, 14, 19, 23, 24) ~ "R2",
  clusters %in% c(7) ~ "R3",
  clusters %in% c(5, 6, 16) ~ "R4",
  clusters %in% c(18) ~ "crop",
  TRUE ~ as.character(clusters))
seurat@meta.data$compartment = compartment
Idents(seurat) = seurat@meta.data$compartment
compartment.averages = AverageExpression(seurat, return.seurat=TRUE, add.ident="sample")
DoHeatmap(compartment.averages, features=compartments, size=2)

p5 = DimPlot(seurat, group.by="compartment", label=TRUE)
p5

It looks like we are able to identify compartments cells are from across both technologies, which is very good.

knitr::kable(seurat@meta.data %>%
  group_by(compartment, technology) %>%
  summarise(n=n()))
compartment technology n
10 10x 76
10 inDrop 323
13 10x 161
13 inDrop 113
17 10x 165
17 inDrop 85
2 10x 100
2 inDrop 683
20 10x 35
20 inDrop 138
21 10x 140
21 inDrop 4
22 10x 10
22 inDrop 77
3 10x 1
3 inDrop 639
8 10x 252
8 inDrop 208
9 10x 106
9 inDrop 328
crop 10x 31
crop inDrop 214
R1 10x 134
R1 inDrop 446
R2 10x 913
R2 inDrop 3403
R3 10x 206
R3 inDrop 266
R4 10x 649
R4 inDrop 699

write seurat

saveRDS(seurat, file.path("results", "seurat-rough.rds"))

integration

Here we wil use the integration functions in Seurat to integrate the 10x and the inDrop datasets.

Idents(seurat) = seurat@meta.data$seurat_clusters
indrop = SubsetData(seurat, subset.name="technology", accept.value="inDrop")
indrop = NormalizeData(object = indrop)
indrop = ScaleData(object = indrop)
indrop = FindVariableFeatures(indrop, do.plot = FALSE)
tenx = SubsetData(seurat, subset.name="technology", accept.value="10x")
tenx = NormalizeData(object = tenx)
tenx = ScaleData(object = tenx)
tenx = FindVariableFeatures(tenx, do.plot = FALSE)
hvf = union(VariableFeatures(indrop), VariableFeatures(tenx))
anchors = FindIntegrationAnchors(object.list=c(indrop, tenx), dims=1:30)
integrated = IntegrateData(anchorset = anchors, dims = 1:30, features.to.integrate=rownames(indrop))
DefaultAssay(integrated) = "integrated"

integrated clustering

Here we use the integrated data to re-do the clustering. After integration, we can see that the 10x and the inDrop datasets overlap with each other much better. We can also get a better idea about what the identities of cells might be.

integrated = ScaleData(integrated, verbose = FALSE)
integrated = RunPCA(integrated, npcs = 30, verbose = FALSE)
integrated = RunUMAP(integrated, reduction = "pca", dims = 1:30)
integrated = FindNeighbors(integrated, reduction = "pca", dims = 1:30)
integrated = FindClusters(integrated, resolution = 0.5)
## Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck
## 
## Number of nodes: 10605
## Number of edges: 407093
## 
## Running Louvain algorithm...
## Maximum modularity in 10 random starts: 0.9284
## Number of communities: 22
## Elapsed time: 0 seconds
p2 = DimPlot(integrated, reduction = "umap", group.by = "technology")
p4 = DimPlot(integrated, reduction = "umap", group.by = "rough_celltype", label=TRUE)
p6 = DimPlot(integrated, group.by="compartment", label=TRUE)
library(gridExtra)
grid.arrange(p1 + ggtitle("before integration"), p2 + ggtitle("after integration"))

grid.arrange(p3 + ggtitle("before integration"), p4 + ggtitle("after integration"))

grid.arrange(p5 + ggtitle("before integration"), p6 + ggtitle("after integration"))

Here we redo the marker expression, using the integrated cluster IDs rather than the unintegrated IDs from before.

DimPlot(integrated, group.by="seurat_clusters", label=TRUE)

integrated.averages = AverageExpression(integrated, return.seurat=TRUE, add.ident="sample")
DoHeatmap(integrated.averages, features=head(markers$genes, 20), size=2)

DoHeatmap(integrated.averages, features=tail(markers$genes, 21), size=2)

knitr::kable(markers)
cell_identity genes
ISC Dl
ISC Smvt
ISC CG10006
pISC sna
ISCR2R5 polo
ISC/EB esg
EB E(spl)m3-HLH
EB E(spl)malpha-BFM
EB E(spl)mbeta-HLH
EB klu
EE pros
EE 7B2
EE_subtype AstA
EE_subtype AstC
EE_subtype CCHa1
EE_subtype CCHa2
EE_subtype Tk
EE_subtype NPF
EE_subtype Orcokinin
EE_subtype Dh31
aEC nub
aEC Myo31DF
mEC nub
mEC Myo31DF
pEC nub
pEC Myo31DF
aEC betaTry
pEC lambdaTry
mEC lab
mEC Vha100-4
cardia Pgant4
copper PGRP-LA
copper PGRP-LC
copper Apoltp
VM Vn
LFC Ilp3
LFC PGRP-SC1a
LFC PGRP-SC1b
LFC PGRP-SD
iron PGRP-SC2
iron ZIP1(Zip42C.1)

Here we are missing classifying cluster 9, 13, and 16. classying 19 and 20 as ISC cells is pretty weak, IMO, so maybe those too.

clusters = integrated@meta.data$seurat_clusters
integrated_celltype = case_when(
  clusters %in% c(3) ~ "ISC/EB",
  clusters %in% c(12, 14, 15) ~ "EE",
  clusters %in% c(1, 4, 18) ~ "aEC",
  clusters %in% c(6) ~ "iron/copper",
  clusters %in% c(5, 11, 13) ~ "pEC",
  clusters %in% c(8) ~ "mEC",
  clusters %in% c(7) ~ "LFC",
  clusters %in% c(17) ~ "cardia",
  TRUE ~ as.character(clusters))
integrated@meta.data$integrated_celltype = integrated_celltype
Idents(integrated) = integrated@meta.data$integrated_celltype
integratedcelltype.averages = AverageExpression(integrated, return.seurat=TRUE, add.ident="sample")
DoHeatmap(integratedcelltype.averages, features=markers$genes, size=2)

p7 = DimPlot(integrated, group.by="integrated_celltype", label=TRUE)
grid.arrange(p4 + ggtitle("integrated celltypes"), p7 + ggtitle("integrated celltypes"))

knitr::kable(integrated@meta.data %>%
             dplyr::select(seurat_clusters, integrated_celltype) %>%
             unique(), row.names=FALSE)
seurat_clusters integrated_celltype
7 LFC
5 pEC
2 2
10 10
0 0
11 pEC
6 iron/copper
19 19
8 mEC
21 21
1 aEC
4 aEC
13 pEC
12 EE
3 ISC/EB
14 EE
9 9
16 16
18 aEC
20 20
15 EE
17 cardia

cluster markers (ROC)

Here we will find markers for each of the clusters where we haven’t been able to assign a cell type so we can hopefully figure out what type of cells they are.

Idents(integrated) = integrated@meta.data$seurat_cluster

Cluster 0

dir.create(file.path("results", "markers_roc"))
markersroc_0 = FindMarkers(integrated, ident.1=0, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markersroc_0, file.path("results", "markers_roc", "cluster0-markers_roc.csv"))
knitr::kable(head(markersroc_0, 25))
gene myAUC avg_diff power pct.1 pct.2
CG12374 0.815 0.7426282 0.630 1.000 0.604
Jon99Ciii 0.790 0.3977139 0.580 0.999 0.669
Bace 0.777 0.4316666 0.554 0.996 0.599
betaTry 0.767 0.5517806 0.534 1.000 0.641
Jon99Cii 0.766 0.3550315 0.532 0.975 0.574
alphaTry 0.759 0.5690638 0.518 1.000 0.790
Jon65Aiii 0.745 0.2787393 0.490 0.975 0.627
CG30025 0.743 0.4160439 0.486 0.973 0.601
pre-rRNA:CR45856 0.257 -1.5008458 0.486 0.133 0.607
Pebp1 0.734 0.3555240 0.468 0.876 0.524
18SrRNA:CR45841 0.273 -1.1015931 0.454 0.158 0.613
18SrRNA-Psi:CR45861 0.273 -1.5898494 0.454 0.091 0.528
18SrRNA:CR41548 0.277 -1.0260107 0.446 0.170 0.623
18SrRNA:CR45838 0.278 -1.0210926 0.444 0.171 0.623
CG12057 0.721 0.6674180 0.442 0.624 0.298
18SrRNA-Psi:CR41602 0.279 -1.5761435 0.442 0.089 0.514
pre-rRNA:CR45846 0.284 -0.9791249 0.432 0.243 0.699
mt:CoI 0.703 0.5018172 0.406 0.984 0.893
mt:CoIII 0.697 0.4695789 0.394 0.961 0.818
His3.3B 0.310 -0.4862776 0.380 0.268 0.682
lncRNA:CR42491 0.311 -0.3863641 0.378 0.189 0.582
14-3-3epsilon 0.316 -0.6482983 0.368 0.169 0.559
crc 0.317 -0.5028142 0.366 0.115 0.479
asRNA:CR45874 0.322 -0.5536836 0.356 0.085 0.421
Hr4 0.323 -0.5593047 0.354 0.142 0.504

Cluster 1

dir.create(file.path("results", "markers_roc"))
markers_roc1 = FindMarkers(integrated, ident.1=1, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc1, file.path("results", "markers_roc", "cluster1-markers_roc.csv"))
knitr::kable(head(markers_roc1, 25))
gene myAUC avg_diff power pct.1 pct.2
CG6295 0.951 2.4180003 0.902 0.970 0.349
CG12374 0.933 1.7278512 0.866 0.991 0.659
betaTry 0.932 1.4927853 0.864 0.992 0.691
CG30025 0.923 1.5042412 0.846 0.992 0.650
alphaTry 0.919 1.2832729 0.838 0.992 0.819
Jon99Ciii 0.893 1.3578745 0.786 0.989 0.715
Jon99Cii 0.877 1.2946068 0.754 0.978 0.629
CG30031 0.869 1.4176099 0.738 0.954 0.505
gammaTry 0.869 1.4176099 0.738 0.954 0.505
deltaTry 0.861 1.3654433 0.722 0.947 0.502
CG5107 0.853 1.4549836 0.706 0.935 0.488
epsilonTry 0.846 1.4368083 0.692 0.934 0.433
Jon65Aiii 0.824 1.2721982 0.648 0.968 0.675
Bace 0.821 0.9766733 0.642 0.979 0.656
Diedel3 0.809 1.2266437 0.618 0.879 0.467
CG17192 0.808 2.2434487 0.616 0.692 0.104
Jon99Ci 0.800 1.3032560 0.600 0.822 0.325
Amy-p 0.797 0.8226741 0.594 0.848 0.395
Jon25Bi 0.796 2.4973360 0.592 0.712 0.171
Jon65Aiv 0.792 1.0718509 0.584 0.966 0.674
tobi 0.780 1.6524216 0.560 0.749 0.284
Mal-A6 0.763 1.6997284 0.526 0.706 0.250
Acbp5 0.242 -2.2991966 0.516 0.396 0.714
rha 0.747 1.2445401 0.494 0.634 0.189
Mal-A1 0.743 0.9628729 0.486 0.827 0.463

Cluster 2

dir.create(file.path("results", "markers_roc"))
markers_roc2 = FindMarkers(integrated, ident.1=2, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc2, file.path("results", "markers_roc", "cluster2-markers_roc.csv"))
knitr::kable(head(markers_roc2, 25))
gene myAUC avg_diff power pct.1 pct.2
CG13315 0.781 1.3570550 0.562 0.888 0.628
pre-rRNA:CR45847 0.234 -1.6052927 0.532 0.435 0.771
MtnC 0.758 0.9374882 0.516 0.911 0.673
pre-rRNA:CR45845 0.248 -1.3662956 0.504 0.546 0.838
RpLP1 0.741 0.6907439 0.482 0.939 0.844
RpLP2 0.741 0.5753640 0.482 0.988 0.948
CG34330 0.738 0.8335108 0.476 0.868 0.654
RpL41 0.738 0.6366086 0.476 0.990 0.972
RpS26 0.734 0.6451404 0.468 0.942 0.857
Cam 0.730 0.9412521 0.460 0.872 0.708
RpL37A 0.728 0.6184733 0.456 0.965 0.914
RpL35 0.727 0.6034513 0.454 0.965 0.889
Acbp5 0.724 1.3897104 0.448 0.870 0.672
RpL19 0.724 0.6111552 0.448 0.950 0.862
RpS21 0.723 0.7083381 0.446 0.885 0.773
RpS30 0.720 0.6082303 0.440 0.935 0.857
CG43349 0.719 1.3817990 0.438 0.643 0.296
RpL23 0.717 0.5709788 0.434 0.965 0.920
sta 0.715 0.5939545 0.430 0.949 0.879
RpL15 0.712 0.6074820 0.424 0.948 0.853
RpS27 0.711 0.5947775 0.422 0.954 0.885
pre-rRNA:CR45846 0.290 -1.8422695 0.420 0.354 0.625
RpS28b 0.710 0.6617437 0.420 0.903 0.800
RpL13 0.707 0.6758232 0.414 0.871 0.746
RpL27A 0.706 0.5697464 0.412 0.933 0.863

Cluster 3

markers_roc3 = FindMarkers(integrated, ident.1=3, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc3, file.path("results", "markers_roc", "cluster3-markers_roc.csv"))
knitr::kable(head(markers_roc3, 25))
gene myAUC avg_diff power pct.1 pct.2
lncRNA:CR40469 0.971 2.3022701 0.942 0.999 0.907
lncRNA:CR34335 0.969 2.0554202 0.938 0.996 0.804
alphaTry 0.114 -3.5902805 0.772 0.302 0.877
Vha16-1 0.143 -2.2895249 0.714 0.243 0.837
betaTry 0.149 -3.7497180 0.702 0.116 0.765
Jon65Aiv 0.151 -3.2747248 0.698 0.089 0.748
Jon99Ciii 0.161 -3.7708719 0.678 0.230 0.779
MtnC 0.168 -2.7958055 0.664 0.108 0.738
Bace 0.169 -3.3875524 0.662 0.118 0.728
CG30025 0.172 -3.6942286 0.656 0.121 0.724
CG13315 0.185 -2.2055407 0.630 0.068 0.695
Jon65Aiii 0.185 -2.9903301 0.630 0.205 0.740
CG12374 0.186 -3.4517539 0.628 0.180 0.727
Acbp5 0.188 -2.5010183 0.624 0.159 0.730
RpL41 0.808 0.5931470 0.616 0.992 0.972
Vha13 0.195 -2.0011931 0.610 0.155 0.724
MtnA 0.201 -2.7763653 0.598 0.260 0.744
CG34330 0.202 -1.9493609 0.596 0.176 0.711
CG44008 0.202 -2.0279360 0.596 0.038 0.626
Jon99Cii 0.203 -3.1741526 0.594 0.155 0.699
cib 0.793 1.8690538 0.586 0.716 0.303
bun 0.793 1.4059197 0.586 0.809 0.419
aqz 0.787 0.9852902 0.574 0.869 0.488
Pebp1 0.223 -2.7887394 0.554 0.095 0.638
Tet 0.774 1.9552201 0.548 0.607 0.167

Cluster 4

markers_roc4 = FindMarkers(integrated, ident.1=4, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc4, file.path("results", "markers_roc", "cluster4-markers_roc.csv"))
knitr::kable(head(markers_roc4, 25))
gene myAUC avg_diff power pct.1 pct.2
Bace 0.948 1.6662209 0.896 1.000 0.657
betaTry 0.915 1.3489190 0.830 1.000 0.693
alphaTry 0.914 1.2351138 0.828 1.000 0.820
CG30025 0.908 1.3850072 0.816 1.000 0.652
CG7542 0.902 1.9216820 0.804 0.917 0.277
Diedel3 0.890 1.4456951 0.780 0.986 0.462
Mal-A7 0.874 2.5066108 0.748 0.846 0.268
CG5107 0.873 1.3861378 0.746 0.976 0.488
epsilonTry 0.865 1.3547922 0.730 0.967 0.434
deltaTry 0.861 1.3768449 0.722 0.974 0.504
CG30031 0.861 1.3478686 0.722 0.972 0.507
gammaTry 0.861 1.3478686 0.722 0.972 0.507
Jon65Aiv 0.851 1.0712409 0.702 0.992 0.674
Mal-A1 0.850 1.6094906 0.700 0.919 0.459
Pebp1 0.843 1.2475968 0.686 0.974 0.567
CG8834 0.833 1.2595427 0.666 0.873 0.261
Jon65Aiii 0.825 0.8930227 0.650 0.991 0.676
CG4377 0.812 1.3900148 0.624 0.820 0.302
Cyp4e1 0.808 1.5399858 0.616 0.788 0.314
CG12374 0.777 0.8567660 0.554 0.946 0.665
tobi 0.774 1.5496700 0.548 0.744 0.289
Mal-A8 0.766 1.4995523 0.532 0.663 0.198
MtnA 0.239 -2.7459049 0.522 0.410 0.731
CG4363 0.756 1.3258361 0.512 0.683 0.231
Mal-A6 0.745 1.1122476 0.490 0.677 0.256

Cluster 5

markers_roc5 = FindMarkers(integrated, ident.1=5, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc5, file.path("results", "markers_roc", "cluster5-markers_roc.csv"))
knitr::kable(head(markers_roc5, 25))
gene myAUC avg_diff power pct.1 pct.2
LManVI 0.977 3.2826160 0.954 0.979 0.243
LManV 0.942 2.9089739 0.884 0.915 0.111
CG31343 0.941 2.0399043 0.882 0.972 0.294
ninaD 0.906 2.6390118 0.812 0.871 0.126
LManII 0.904 2.1743168 0.808 0.869 0.141
CG15534 0.887 2.2950064 0.774 0.852 0.135
iotaTry 0.868 1.8477332 0.736 0.836 0.170
Cyp6d5 0.865 1.6516791 0.730 0.891 0.332
Mal-A4 0.859 1.9131132 0.718 0.813 0.161
CG11911 0.855 1.2775764 0.710 0.940 0.346
lambdaTry 0.851 1.5033115 0.702 0.862 0.217
Ance 0.836 1.2035410 0.672 0.882 0.282
LManI 0.834 1.7075755 0.668 0.760 0.109
CG7025 0.834 1.7006117 0.668 0.792 0.132
CG14629 0.831 1.4508569 0.662 0.741 0.135
asRNA:CR45281 0.825 1.1363389 0.650 0.746 0.158
CG9673 0.823 1.0930344 0.646 0.848 0.220
LManIII 0.821 2.2069193 0.642 0.697 0.052
CG15533 0.814 1.4710121 0.628 0.688 0.063
fabp 0.806 0.9226544 0.612 0.929 0.461
CG33127 0.800 1.1443035 0.600 0.836 0.275
Gba1a 0.800 1.1219500 0.600 0.660 0.082
CG11151 0.800 0.9555007 0.600 0.859 0.307
CHMP2B 0.800 0.9366572 0.600 0.704 0.218
CG33306 0.798 1.1217776 0.596 0.718 0.123

Cluster 6

markers_roc6 = FindMarkers(integrated, ident.1=6, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc6, file.path("results", "markers_roc", "cluster6-markers_roc.csv"))
knitr::kable(head(markers_roc6, 25))
gene myAUC avg_diff power pct.1 pct.2
MtnC 0.882 1.7872526 0.764 0.973 0.675
Arc1 0.867 2.4032135 0.734 0.911 0.374
CG43774 0.855 3.1529264 0.710 0.807 0.207
CG34330 0.838 1.4639846 0.676 0.945 0.655
CG5399 0.837 1.6595736 0.674 0.857 0.351
mbl 0.834 2.1404858 0.668 0.850 0.388
MtnD 0.812 2.6718753 0.624 0.775 0.307
MtnA 0.811 1.4680532 0.622 0.959 0.693
CG15423 0.807 1.9412040 0.614 0.737 0.166
CG15422 0.776 1.7420240 0.552 0.730 0.245
Stat92E 0.769 1.7067355 0.538 0.631 0.212
Vha16-1 0.754 0.9391610 0.508 0.952 0.783
MtnB 0.750 1.6650008 0.500 0.708 0.303
Vha13 0.739 0.9406012 0.478 0.900 0.669
Adat1 0.738 1.0285194 0.476 0.640 0.249
thetaTry 0.730 1.9437544 0.460 0.623 0.136
CycG 0.726 1.2622215 0.452 0.812 0.503
MtnE 0.723 1.6779141 0.446 0.642 0.249
pre-rRNA:CR45845 0.279 -1.2473915 0.442 0.621 0.826
RpLP2 0.720 0.6403127 0.440 0.986 0.949
pre-rRNA:CR45847 0.281 -1.3870317 0.438 0.571 0.755
CG8177 0.715 2.2074678 0.430 0.574 0.251
Msp300 0.709 1.2365695 0.418 0.719 0.365
Lk6 0.705 1.2882895 0.410 0.769 0.492
RpL41 0.705 0.5436570 0.410 0.998 0.972

Cluster 7

markers_roc7 = FindMarkers(integrated, ident.1=7, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc7, file.path("results", "markers_roc", "cluster7-markers_roc.csv"))
knitr::kable(head(markers_roc7, 25))
gene myAUC avg_diff power pct.1 pct.2
CG10472 0.984 2.616301 0.968 1.000 0.405
Jon99Fii 0.934 3.236508 0.868 0.940 0.255
Jon25Biii 0.926 3.256755 0.852 0.960 0.411
CG17571 0.922 2.704814 0.844 0.936 0.257
mag 0.922 2.295322 0.844 0.929 0.172
Acbp5 0.917 1.379220 0.834 0.998 0.670
CG11911 0.904 2.168190 0.808 0.933 0.348
Jon65Ai 0.896 2.885554 0.792 0.858 0.142
Jon99Fi 0.880 3.266696 0.760 0.862 0.228
CG7953 0.861 2.720265 0.722 0.800 0.203
CG16749 0.858 1.553809 0.716 0.894 0.276
CG17633 0.855 2.406228 0.710 0.811 0.222
CG3868 0.849 1.400003 0.698 0.902 0.330
CG8997 0.845 3.178884 0.690 0.801 0.263
Jon65Aii 0.836 3.293820 0.672 0.814 0.210
lncRNA:CR40469 0.170 -2.088592 0.660 0.718 0.925
CG31198 0.829 1.287501 0.658 0.885 0.290
CG7916 0.821 3.257496 0.642 0.750 0.252
Jon99Cii 0.821 1.725367 0.642 0.954 0.642
CG18493 0.820 1.598612 0.640 0.794 0.239
CG31323 0.810 1.136244 0.620 0.836 0.322
Jon99Ciii 0.809 1.736809 0.618 0.965 0.725
CG8093 0.802 1.531794 0.604 0.670 0.092
asRNA:CR45281 0.802 1.095471 0.604 0.709 0.161
CG31233 0.799 1.109342 0.598 0.829 0.263

Cluster 8

markers_roc8 = FindMarkers(integrated, ident.1=8, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc8, file.path("results", "markers_roc", "cluster8-markers_roc.csv"))
knitr::kable(head(markers_roc8, 25))
gene myAUC avg_diff power pct.1 pct.2
Vha16-1 0.965 2.334783 0.930 1.000 0.783
CG5767 0.963 4.182364 0.926 0.949 0.170
Vha13 0.942 1.866896 0.884 0.996 0.667
MtnB 0.940 3.424423 0.880 0.947 0.297
CG30479 0.914 3.367492 0.828 0.898 0.135
Adat1 0.909 2.418973 0.818 0.889 0.242
CG15423 0.907 2.526331 0.814 0.895 0.165
Vha55 0.899 2.306491 0.798 0.958 0.416
CG30480 0.893 2.297150 0.786 0.802 0.056
Vha68-2 0.892 2.124161 0.784 0.931 0.341
CG5399 0.887 2.160550 0.774 0.913 0.354
Vha44 0.886 2.244559 0.772 0.964 0.474
Vha14-1 0.883 1.918986 0.766 0.900 0.299
VhaAC39-1 0.881 2.244777 0.762 0.913 0.328
CG15422 0.868 1.932872 0.736 0.884 0.243
Tsp42Ec 0.865 2.502074 0.730 0.833 0.179
VhaSFD 0.864 2.068645 0.728 0.869 0.244
Vha26 0.849 1.874010 0.698 0.878 0.308
Vha36-1 0.845 1.801660 0.690 0.862 0.304
Argk 0.844 1.768606 0.688 0.909 0.389
VhaAC45 0.838 1.735004 0.676 0.878 0.320
Mpcp2 0.833 1.875708 0.666 0.951 0.558
VhaM9.7-b 0.832 1.618574 0.664 0.842 0.307
CG7430 0.827 2.010335 0.654 0.793 0.292
CAH1 0.824 2.278490 0.648 0.777 0.132

Cluster 9

markers_roc9 = FindMarkers(integrated, ident.1=9, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc9, file.path("results", "markers_roc", "cluster9-markers_roc.csv"))
knitr::kable(head(markers_roc9, 25))
gene myAUC avg_diff power pct.1 pct.2
Amy-p 0.927 3.7546156 0.854 0.926 0.413
Amy-d 0.823 3.4545206 0.646 0.718 0.233
CG11400 0.820 3.3465943 0.640 0.688 0.146
lncRNA:CR40469 0.811 0.6702035 0.622 0.995 0.911
CG4928 0.804 3.3092892 0.608 0.678 0.197
CG3819 0.803 3.3652770 0.606 0.657 0.126
CG15043 0.793 3.0710478 0.586 0.642 0.159
lncRNA:CR34335 0.781 0.7487954 0.562 0.980 0.813
CG6839 0.766 3.6581995 0.532 0.622 0.156
trol 0.755 0.7098187 0.510 0.995 0.990
fmt 0.749 0.6997454 0.498 0.997 0.995
CG14125 0.742 4.9471279 0.484 0.563 0.197
CG13323 0.740 2.4128346 0.480 0.612 0.282
Mal-A1 0.283 -1.8074388 0.434 0.124 0.508
Vha13 0.712 0.6736012 0.424 0.883 0.673
28SrRNA-Psi:CR40741 0.292 -1.5479371 0.416 0.124 0.464
Adhr 0.708 1.6133663 0.416 0.541 0.202
alphaTry 0.295 -1.7875947 0.410 0.761 0.836
CG30025 0.295 -1.8352758 0.410 0.393 0.689
CG15818 0.704 3.2587424 0.408 0.492 0.108
sesB 0.703 0.6889130 0.406 0.893 0.772
Oda 0.699 0.7258382 0.398 0.919 0.830
pAbp 0.695 0.5943689 0.390 0.898 0.761
Vha16-1 0.692 0.2603364 0.384 0.942 0.786
Adh 0.690 1.5147031 0.380 0.503 0.189

Cluster 10

markers_roc10 = FindMarkers(integrated, ident.1=10, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc10, file.path("results", "markers_roc", "cluster10-markers_roc.csv"))
knitr::kable(head(markers_roc10, 25))
gene myAUC avg_diff power pct.1 pct.2
pre-rRNA:CR45846 0.883 1.8543423 0.766 0.967 0.592
pre-rRNA:CR45847 0.877 1.6122745 0.754 0.994 0.737
18SrRNA:CR41548 0.875 1.8092135 0.750 0.917 0.516
18SrRNA:CR45838 0.875 1.8088839 0.750 0.917 0.516
18SrRNA:CR45841 0.869 1.8323186 0.738 0.902 0.506
pre-rRNA:CR45845 0.833 1.5562434 0.666 0.979 0.810
pre-rRNA:CR45856 0.828 1.9342579 0.656 0.842 0.497
18SrRNA-Psi:CR45861 0.807 1.9929792 0.614 0.774 0.425
18SrRNA-Psi:CR41602 0.768 1.9522780 0.536 0.714 0.415
28SrRNA:CR45844 0.715 1.1993126 0.430 0.744 0.607
28SrRNA-Psi:CR45851 0.702 1.3256885 0.404 0.688 0.531
Jon99Cii 0.702 0.5646288 0.404 0.893 0.650
Jon99Ciii 0.681 0.2742590 0.362 0.938 0.731
CG12374 0.678 0.3401385 0.356 0.943 0.678
Bace 0.670 0.2581049 0.340 0.923 0.674
CG30025 0.662 0.2913271 0.324 0.899 0.671
14-3-3epsilon 0.347 -0.6633088 0.306 0.158 0.488
28SrRNA-Psi:CR40741 0.649 1.3577011 0.298 0.565 0.448
28SrRNA-Psi:CR41609 0.642 1.3964064 0.284 0.560 0.459
clos 0.364 -0.3143394 0.272 0.107 0.384
CG33995 0.366 -0.3927518 0.268 0.071 0.341
crc 0.366 -0.3928422 0.268 0.131 0.412
His3.3B 0.366 -0.4291378 0.268 0.304 0.605
CG10472 0.633 0.3546737 0.266 0.601 0.430
CG42322 0.369 -0.2620372 0.262 0.182 0.471

Cluster 11

markers_roc11 = FindMarkers(integrated, ident.1=11, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc11, file.path("results", "markers_roc", "cluster11-markers_roc.csv"))
knitr::kable(head(markers_roc11, 25))
gene myAUC avg_diff power pct.1 pct.2
Gs2 0.960 2.5571275 0.920 0.978 0.369
Acbp3 0.957 2.7989831 0.914 0.978 0.298
CG31343 0.927 1.8074765 0.854 0.984 0.310
zetaTry 0.902 2.2494411 0.804 0.873 0.180
CG13492 0.891 1.8982504 0.782 0.876 0.217
CG9673 0.880 1.7049865 0.760 0.895 0.234
fabp 0.864 1.2838541 0.728 0.962 0.471
CG15254 0.842 2.5365965 0.684 0.796 0.138
MtnA 0.842 1.1704878 0.684 0.987 0.699
CG10116 0.837 1.7645785 0.674 0.812 0.219
Gdh 0.837 1.7188263 0.674 0.768 0.192
CG8774 0.834 2.0324823 0.668 0.768 0.148
Agpat4 0.819 1.5228590 0.638 0.761 0.194
CG32473 0.817 2.6193899 0.634 0.739 0.143
CG5958 0.801 1.7785181 0.602 0.726 0.177
Pepck2 0.795 1.7380777 0.590 0.653 0.096
NAAT1 0.790 1.9975076 0.580 0.669 0.058
dmGlut 0.790 1.4212022 0.580 0.653 0.092
CG10911 0.789 0.9464320 0.578 0.898 0.402
CG32633 0.788 1.2566876 0.576 0.758 0.236
Jon99Ciii 0.215 -3.9366132 0.570 0.455 0.746
Jon65Aiv 0.218 -3.0364897 0.564 0.366 0.708
CG4653 0.781 1.4001194 0.562 0.736 0.220
Acbp5 0.779 0.7779398 0.558 0.978 0.679
Jon65Aiii 0.221 -3.0295947 0.558 0.360 0.710

Cluster 12

markers_roc12 = FindMarkers(integrated, ident.1=12, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc12, file.path("results", "markers_roc", "cluster12-markers_roc.csv"))
knitr::kable(head(markers_roc12, 25))
gene myAUC avg_diff power pct.1 pct.2
IA-2 0.964 2.647187 0.928 0.993 0.403
nrv3 0.916 2.247964 0.832 0.909 0.242
7B2 0.913 2.096422 0.826 0.915 0.231
Phm 0.882 1.793405 0.764 0.866 0.234
Pal2 0.852 2.079188 0.704 0.782 0.147
pros 0.844 2.015971 0.688 0.769 0.181
w 0.841 2.291362 0.682 0.759 0.206
CG30183 0.838 1.964119 0.676 0.736 0.130
alphaTry 0.172 -3.018572 0.656 0.440 0.845
Jon99Ciii 0.179 -3.040954 0.642 0.212 0.753
Hsp83 0.813 1.332282 0.626 0.850 0.374
Jon65Aiv 0.195 -2.829406 0.610 0.186 0.714
betaTry 0.195 -3.142426 0.610 0.215 0.731
heph 0.804 1.849458 0.608 0.664 0.104
Galphao 0.803 1.910813 0.606 0.691 0.175
Jon65Aiii 0.200 -2.671267 0.600 0.199 0.715
RpL8 0.201 -1.069680 0.598 0.560 0.903
Gbeta13F 0.798 1.369896 0.596 0.779 0.335
Bace 0.202 -2.942945 0.596 0.160 0.698
RpL35 0.203 -1.058824 0.594 0.564 0.904
RpL37A 0.207 -1.029812 0.586 0.664 0.925
CG12374 0.207 -3.043524 0.586 0.189 0.701
RpL17 0.211 -1.002210 0.578 0.560 0.895
Jon99Cii 0.216 -2.491872 0.568 0.147 0.673
CG30025 0.216 -3.165431 0.568 0.238 0.691

Cluster 13

markers_roc13 = FindMarkers(integrated, ident.1=13, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc13, file.path("results", "markers_roc", "cluster13-markers_roc.csv"))
knitr::kable(head(markers_roc13, 25))
gene myAUC avg_diff power pct.1 pct.2
CG33926 0.924 2.1060157 0.848 0.950 0.228
CG10912 0.911 2.2182013 0.822 0.907 0.203
CG10911 0.891 1.5862509 0.782 0.963 0.401
CG32368 0.883 2.1692952 0.766 0.854 0.168
Mur29B 0.855 2.5662759 0.710 0.804 0.184
CG42825 0.844 2.2269327 0.688 0.761 0.116
CG14499 0.843 3.0172332 0.686 0.764 0.075
lectin-37Da 0.839 2.7082242 0.678 0.757 0.081
CG31086 0.823 1.5386260 0.646 0.804 0.311
CG31323 0.821 1.7430702 0.642 0.837 0.335
CG18493 0.794 1.4654841 0.588 0.754 0.253
RpL3 0.213 -0.9937845 0.574 0.718 0.896
eEF1alpha1 0.218 -0.8581779 0.564 0.887 0.970
Npc2e 0.778 2.1977119 0.556 0.691 0.154
CG15255 0.778 1.3555371 0.556 0.771 0.267
asRNA:CR44192 0.778 1.3426450 0.556 0.635 0.150
Fer2LCH 0.772 0.9810805 0.544 0.897 0.479
RpS29 0.238 -0.8599220 0.524 0.837 0.941
lectin-37Db 0.757 1.3335771 0.514 0.688 0.226
CG9568 0.754 1.6351462 0.508 0.688 0.224
CG44008 0.754 0.6956464 0.508 0.924 0.571
RpS2 0.247 -0.9726861 0.506 0.731 0.873
sta 0.249 -0.9397450 0.502 0.721 0.889
RpL14 0.250 -0.9298185 0.500 0.638 0.842
RpL38 0.254 -0.9220046 0.492 0.585 0.843

Cluster 14

markers_roc14 = FindMarkers(integrated, ident.1=14, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc14, file.path("results", "markers_roc", "cluster14-markers_roc.csv"))
knitr::kable(head(markers_roc14, 25))
gene myAUC avg_diff power pct.1 pct.2
NPF 0.989 5.147902 0.978 0.983 0.290
IA-2 0.986 3.213505 0.972 1.000 0.403
7B2 0.933 3.044842 0.866 0.924 0.231
Phm 0.913 2.665425 0.826 0.890 0.234
Tk 0.898 2.448163 0.796 0.857 0.210
svr 0.896 2.567417 0.792 0.867 0.243
chrb 0.884 3.500329 0.768 0.841 0.281
nrv3 0.867 2.251026 0.734 0.827 0.245
CG30183 0.852 2.181997 0.704 0.764 0.129
CG46385 0.850 2.210807 0.700 0.831 0.325
esg 0.840 2.222163 0.680 0.761 0.184
unc-13-4A 0.834 1.715761 0.668 0.757 0.169
Hsp23 0.827 2.841930 0.654 0.781 0.318
Ldh 0.827 2.230674 0.654 0.728 0.177
alphaTry 0.177 -2.887524 0.646 0.482 0.844
cpo 0.822 2.446349 0.644 0.761 0.275
Pal2 0.822 2.156723 0.644 0.748 0.149
Hsp22 0.822 1.852289 0.644 0.777 0.341
mbl 0.821 1.258315 0.642 0.867 0.399
Hsp83 0.820 1.744655 0.640 0.847 0.374
Unr 0.815 1.970741 0.630 0.804 0.377
Hsp26 0.809 1.839010 0.618 0.847 0.432
Ih 0.807 2.434480 0.614 0.661 0.097
shep 0.807 1.896214 0.614 0.728 0.241
CG30025 0.195 -2.890765 0.610 0.199 0.692

Cluster 15

markers_roc15 = FindMarkers(integrated, ident.1=15, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc15, file.path("results", "markers_roc", "cluster15-markers_roc.csv"))
knitr::kable(head(markers_roc15, 25))
gene myAUC avg_diff power pct.1 pct.2
Orcokinin 1.000 5.327171 1.000 1.000 0.292
AstC 0.987 3.333166 0.974 0.996 0.293
svr 0.965 3.009325 0.930 0.979 0.244
unc-13-4A 0.957 2.946501 0.914 0.963 0.167
nrv3 0.939 2.642404 0.878 0.954 0.246
IA-2 0.936 2.182955 0.872 0.992 0.407
7B2 0.927 2.126905 0.854 0.959 0.235
Phm 0.923 2.484770 0.846 0.934 0.237
tap 0.917 1.980146 0.834 0.851 0.030
heph 0.898 2.319493 0.796 0.846 0.103
CG30183 0.898 1.741654 0.796 0.880 0.130
pros 0.896 2.421778 0.792 0.871 0.183
CG14989 0.884 3.799407 0.768 0.871 0.190
Ca-alpha1T 0.880 2.284545 0.760 0.788 0.055
slo 0.880 1.877297 0.760 0.780 0.045
Rgk3 0.879 2.170133 0.758 0.768 0.022
jus 0.873 1.609255 0.746 0.768 0.051
Pal2 0.871 2.284061 0.742 0.846 0.150
epsilonTry 0.129 -2.541908 0.742 0.046 0.484
lncRNA:CR31451 0.870 1.543013 0.740 0.780 0.062
fkh 0.867 1.852307 0.734 0.871 0.279
CG44247 0.867 1.801113 0.734 0.780 0.060
shep 0.864 1.872548 0.728 0.876 0.241
CG46385 0.859 1.946008 0.718 0.905 0.327
fru 0.859 1.606407 0.718 0.780 0.144

Cluster 16

markers_roc16 = FindMarkers(integrated, ident.1=16, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc16, file.path("results", "markers_roc", "cluster16-markers_roc.csv"))
knitr::kable(head(markers_roc16, 25))
gene myAUC avg_diff power pct.1 pct.2
EbpIII 0.948 4.834704 0.896 0.920 0.223
Acbp2 0.914 2.549511 0.828 0.903 0.364
lncRNA:CR34335 0.905 1.580156 0.810 1.000 0.815
Adhr 0.883 3.095023 0.766 0.819 0.200
Adh 0.861 3.083127 0.722 0.776 0.187
lncRNA:CR40469 0.858 1.202531 0.716 1.000 0.912
to 0.853 3.587573 0.706 0.726 0.064
spidey 0.841 3.329530 0.682 0.734 0.216
CG30197 0.825 3.532219 0.650 0.684 0.094
CG31523 0.824 2.969051 0.648 0.658 0.148
CG31522 0.807 3.284502 0.614 0.688 0.226
CG1124 0.786 3.623984 0.572 0.620 0.110
Msr-110 0.779 2.666640 0.558 0.612 0.116
CG6770 0.768 1.183859 0.536 0.819 0.557
ple 0.763 2.881654 0.526 0.553 0.030
wat 0.756 3.231272 0.512 0.532 0.046
CG10237 0.750 2.613266 0.500 0.536 0.048
Moe 0.748 1.735258 0.496 0.591 0.240
emp 0.744 2.438345 0.488 0.523 0.089
vir-1 0.738 2.554232 0.476 0.494 0.033
Hacd1 0.737 2.129050 0.474 0.485 0.111
CG8306 0.734 2.509803 0.468 0.494 0.060
ATPCL 0.730 2.223677 0.460 0.527 0.170
Inos 0.729 2.157946 0.458 0.502 0.099
Taldo 0.727 1.879034 0.454 0.523 0.197

Cluster 17

markers_roc17 = FindMarkers(integrated, ident.1=17, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc17, file.path("results", "markers_roc", "cluster17-markers_roc.csv"))
knitr::kable(head(markers_roc17, 25))
gene myAUC avg_diff power pct.1 pct.2
CG14645 0.999 5.166353 0.998 1.000 0.282
Skp2 0.998 5.300042 0.996 1.000 0.389
CG3906 0.992 5.087155 0.984 0.991 0.217
CG34324 0.990 5.281636 0.980 0.987 0.277
CG34220 0.988 5.093672 0.976 0.987 0.420
CG4783 0.971 3.204050 0.942 0.956 0.160
Muc68D 0.958 5.030528 0.916 0.933 0.207
CG11672 0.932 3.842805 0.864 0.871 0.037
Idgf4 0.874 2.455162 0.748 0.782 0.072
Pgant4 0.828 2.450677 0.656 0.662 0.014
alphaTry 0.186 -2.869711 0.628 0.440 0.842
CG30025 0.186 -3.163239 0.628 0.116 0.690
Jon99Ciii 0.201 -2.909744 0.598 0.227 0.749
CG34330 0.214 -2.012285 0.572 0.133 0.682
Jon65Aiv 0.216 -2.543440 0.568 0.204 0.709
Jon65Aiii 0.217 -2.454261 0.566 0.209 0.710
betaTry 0.217 -2.727593 0.566 0.236 0.726
Jon99Cii 0.229 -2.470175 0.542 0.151 0.669
Bace 0.232 -2.477221 0.536 0.209 0.693
Vha16-1 0.239 -1.611035 0.522 0.364 0.801
MtnC 0.239 -2.133976 0.522 0.280 0.700
sesB 0.240 -1.164468 0.520 0.320 0.787
Act5C 0.243 -1.032861 0.514 0.604 0.900
CG12374 0.245 -2.491602 0.510 0.249 0.696
CG30026 0.754 1.759105 0.508 0.511 0.003

Cluster 18

markers_roc18 = FindMarkers(integrated, ident.1=18, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc18, file.path("results", "markers_roc", "cluster18-markers_roc.csv"))
knitr::kable(head(markers_roc18, 25))
gene myAUC avg_diff power pct.1 pct.2
Npc2f 0.986 3.6101816 0.972 0.988 0.150
Jon65Aiv 0.986 2.1546865 0.972 1.000 0.694
thetaTry 0.980 3.0983381 0.960 0.988 0.149
yip7 0.979 2.3389780 0.958 1.000 0.539
Jon65Aiii 0.979 2.0314714 0.958 1.000 0.695
Peritrophin-15a 0.938 2.6028515 0.876 0.950 0.184
Bace 0.926 1.8774625 0.852 0.994 0.678
CG30031 0.915 1.3828164 0.830 1.000 0.535
gammaTry 0.915 1.3828164 0.830 1.000 0.535
deltaTry 0.890 1.3582514 0.780 1.000 0.532
CG4734 0.885 1.8411515 0.770 0.826 0.132
CG18404 0.879 4.0870292 0.758 0.845 0.161
alphaTry 0.872 1.0091274 0.744 1.000 0.831
Pebp1 0.867 1.5259964 0.734 0.944 0.592
CG9682 0.851 2.5655588 0.702 0.727 0.027
CG30025 0.842 1.0353096 0.684 1.000 0.673
LysB 0.835 1.3202083 0.670 0.826 0.197
CG4830 0.827 2.6994819 0.654 0.789 0.017
lncRNA:CR40469 0.181 -2.0822384 0.638 0.826 0.915
CG4563 0.816 2.4116003 0.632 0.770 0.061
CG5107 0.799 0.8569231 0.598 0.963 0.518
betaTry 0.788 0.8541648 0.576 1.000 0.711
CG34040 0.783 0.5422801 0.566 0.640 0.082
CG4377 0.779 0.9930001 0.558 0.857 0.333
epsilonTry 0.778 0.7924996 0.556 0.957 0.467

Cluster 19

markers_roc19 = FindMarkers(integrated, ident.1=19, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc19, file.path("results", "markers_roc", "cluster19-markers_roc.csv"))
knitr::kable(head(markers_roc19, 25))
gene myAUC avg_diff power pct.1 pct.2
CG6409 0.980 4.8036985 0.960 0.966 0.145
CG42656 0.951 4.4368015 0.902 0.910 0.082
IM33 0.908 3.9650892 0.816 0.831 0.081
Vha13 0.906 1.5702548 0.812 0.978 0.679
lncRNA:CR40469 0.902 1.3818689 0.804 1.000 0.913
Vha16-1 0.892 1.1695161 0.784 0.978 0.791
sesB 0.880 1.2228950 0.760 0.978 0.775
Vha55 0.845 1.5750791 0.690 0.876 0.436
lncRNA:CR34335 0.842 1.0495865 0.684 0.989 0.817
ATPsynC 0.817 1.1686827 0.634 0.910 0.659
Vha44 0.815 1.3559872 0.630 0.843 0.492
Gdh 0.811 2.5818285 0.622 0.685 0.205
Vha26 0.808 1.6921151 0.616 0.764 0.328
Vha68-2 0.805 1.5239685 0.610 0.775 0.362
porin 0.799 1.5868938 0.598 0.809 0.503
tws 0.797 2.6627324 0.594 0.640 0.170
blw 0.790 1.4953967 0.580 0.820 0.556
Ggamma30A 0.778 2.4257360 0.556 0.584 0.061
CG1143 0.776 2.8464808 0.552 0.584 0.052
Irk2 0.775 2.8819032 0.550 0.584 0.057
Msr-110 0.764 2.2743860 0.528 0.596 0.123
pAbp 0.760 0.7985467 0.520 0.910 0.765
Argk 0.753 1.2657596 0.506 0.730 0.408
Mpcp2 0.753 1.1475929 0.506 0.798 0.572
Dic1 0.752 1.8342125 0.504 0.562 0.178

Cluster 20

markers_roc20 = FindMarkers(integrated, ident.1=20, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc20, file.path("results", "markers_roc", "cluster20-markers_roc.csv"))
knitr::kable(head(markers_roc20, 25))
gene myAUC avg_diff power pct.1 pct.2
CG13285 0.980 3.8504561 0.960 0.963 0.020
Msr-110 0.970 3.1788920 0.940 0.963 0.123
CG44013 0.887 2.9005776 0.774 0.778 0.020
lncRNA:CR34335 0.886 1.2155838 0.772 1.000 0.818
Idgf4 0.885 2.3408393 0.770 0.815 0.083
Idgf6 0.881 2.6238014 0.762 0.778 0.043
lncRNA:CR40469 0.876 0.9434807 0.752 1.000 0.913
CG5065 0.852 2.2587751 0.704 0.741 0.095
CG43394 0.837 2.6414936 0.674 0.685 0.028
CG10096 0.829 2.8104611 0.658 0.667 0.012
mt:lrRNA 0.819 0.3111689 0.638 1.000 1.000
Vha16-1 0.186 -2.1373285 0.628 0.241 0.795
MtnC 0.201 -2.6273170 0.598 0.111 0.694
Act5C 0.202 -1.2459382 0.596 0.574 0.896
Orcokinin 0.794 0.5212523 0.588 0.741 0.306
alphaTry 0.208 -2.5242600 0.584 0.444 0.836
Jon65Aiii 0.212 -2.6360940 0.576 0.148 0.702
Jon99Cii 0.227 -2.9229030 0.546 0.093 0.661
GstD11 0.768 2.0589475 0.536 0.537 0.001
Jon65Aiv 0.232 -2.4779429 0.536 0.204 0.701
Acbp5 0.234 -1.7272450 0.532 0.148 0.690
Cys 0.765 1.8452262 0.530 0.611 0.177
CG30025 0.240 -2.5843813 0.520 0.204 0.681
Vha13 0.242 -1.7401562 0.516 0.167 0.684
Atpalpha 0.754 1.4416261 0.508 0.667 0.310

Cluster 21

markers_roc21 = FindMarkers(integrated, ident.1=21, test.use="roc") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(-power)
readr::write_csv(markers_roc21, file.path("results", "markers_roc", "cluster21-markers_roc.csv"))
knitr::kable(head(markers_roc21, 25))
gene myAUC avg_diff power pct.1 pct.2
MtnA 0.933 1.8253924 0.866 1.000 0.706
Alp4 0.810 4.0091710 0.620 0.639 0.055
Vha13 0.801 1.0715504 0.602 0.889 0.681
CG6726 0.784 3.1950727 0.568 0.611 0.109
Mur18B 0.773 3.7634567 0.546 0.583 0.080
Alp2 0.771 3.2290870 0.542 0.556 0.029
sesB 0.756 0.7051319 0.512 0.944 0.776
CG3168 0.746 3.7752656 0.492 0.528 0.073
Vha26 0.745 1.0301148 0.490 0.694 0.330
alphaTry 0.264 -1.4814880 0.472 0.556 0.835
HDAC6 0.729 2.5162694 0.458 0.500 0.161
Mal-A1 0.272 -2.6635592 0.456 0.028 0.495
CG18324 0.723 1.7456577 0.446 0.444 0.126
UGP 0.721 2.2285553 0.442 0.500 0.160
Vha16-1 0.721 0.6729296 0.442 0.889 0.792
CG30025 0.286 -1.8776482 0.428 0.333 0.679
smp-30 0.707 2.1958556 0.414 0.472 0.166
Vha44 0.703 0.8803877 0.406 0.694 0.494
CG31288 0.700 2.0329887 0.400 0.444 0.172
CG30031 0.301 -1.3683171 0.398 0.167 0.543
gammaTry 0.301 -1.3683171 0.398 0.167 0.543
CG13868 0.698 1.8833634 0.396 0.417 0.137
deltaTry 0.303 -1.3374017 0.394 0.167 0.540
CG5107 0.303 -2.0925566 0.394 0.139 0.526
Vha55 0.695 0.8688583 0.390 0.667 0.439

Specific markers-roc plot for cluster 21.

DimPlot(integrated, label=TRUE)

FeaturePlot(integrated, features=c("MtnA", "Vha13", "Smvt", "Oatp58Dc"))

cells21 = rownames(subset(integrated@meta.data, seurat_clusters == 21))
FeaturePlot(integrated[,cells21], features=c("MtnA", "Vha13", "Smvt", "Oatp58Dc")) +
  ggtitle("cluster 21 only")

cluster markers (MAST)

Here we will find markers for each of the clusters where we haven’t been able to assign a cell type so we can hopefully figure out what type of cells they are.

Cluster 0

dir.create(file.path("results", "markers_MAST"))
markers_MAST0 = FindMarkers(integrated, ident.1=0, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST0, file.path("results", "markers_MAST", "cluster0-markers_MAST.csv"))
knitr::kable(head(markers_MAST0, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
CG12374 0 0.7426282 1.000 0.604 0
alphaTry 0 0.5690638 1.000 0.790 0
betaTry 0 0.5517806 1.000 0.641 0
Bace 0 0.4316666 0.996 0.599 0
CG30025 0 0.4160439 0.973 0.601 0
Jon99Ciii 0 0.3977139 0.999 0.669 0
Jon99Cii 0 0.3550315 0.975 0.574 0
Jon65Aiii 0 0.2787393 0.975 0.627 0
pre-rRNA:CR45846 0 -0.9791249 0.243 0.699 0
18SrRNA:CR45838 0 -1.0210926 0.171 0.623 0
18SrRNA:CR41548 0 -1.0260107 0.170 0.623 0
18SrRNA:CR45841 0 -1.1015931 0.158 0.613 0
pre-rRNA:CR45856 0 -1.5008458 0.133 0.607 0
18SrRNA-Psi:CR41602 0 -1.5761435 0.089 0.514 0
18SrRNA-Psi:CR45861 0 -1.5898494 0.091 0.528 0
CG12057 0 0.6674180 0.624 0.298 0
His3.3B 0 -0.4862776 0.268 0.682 0
Pebp1 0 0.3555240 0.876 0.524 0
ps 0 -0.3731687 0.233 0.625 0
lncRNA:CR42491 0 -0.3863641 0.189 0.582 0
14-3-3epsilon 0 -0.6482983 0.169 0.559 0
cib 0 -0.5655328 0.066 0.405 0
CycG 0 -0.8424565 0.207 0.601 0
porin 0 -0.3450815 0.217 0.582 0
28SrRNA-Psi:CR45862 0 -0.2617602 0.030 0.339 0

Cluster 1

dir.create(file.path("results", "markers_MAST"))
markers_MAST1 = FindMarkers(integrated, ident.1=1, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST1, file.path("results", "markers_MAST", "cluster1-markers_MAST.csv"))
knitr::kable(head(markers_MAST1, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
Jon25Bi 0 2.4973360 0.712 0.171 0
CG6295 0 2.4180003 0.970 0.349 0
CG17192 0 2.2434487 0.692 0.104 0
CG34026 0 1.8928975 0.535 0.076 0
CG12374 0 1.7278512 0.991 0.659 0
CG30025 0 1.5042412 0.992 0.650 0
betaTry 0 1.4927853 0.992 0.691 0
CG30031 0 1.4176099 0.954 0.505 0
gammaTry 0 1.4176099 0.954 0.505 0
deltaTry 0 1.3654433 0.947 0.502 0
CG5107 0 1.4549836 0.935 0.488 0
Jon99Cii 0 1.2946068 0.978 0.629 0
epsilonTry 0 1.4368083 0.934 0.433 0
Jon99Ciii 0 1.3578745 0.989 0.715 0
Jon25Bii 0 2.2726902 0.610 0.146 0
alphaTry 0 1.2832729 0.992 0.819 0
Jon99Ci 0 1.3032560 0.822 0.325 0
Mal-A6 0 1.6997284 0.706 0.250 0
Jon65Aiii 0 1.2721982 0.968 0.675 0
Bace 0 0.9766733 0.979 0.656 0
Diedel3 0 1.2266437 0.879 0.467 0
MtnA 0 -2.4875220 0.519 0.724 0
Amy-p 0 0.8226741 0.848 0.395 0
tobi 0 1.6524216 0.749 0.284 0
Acbp5 0 -2.2991966 0.396 0.714 0

Cluster 2

dir.create(file.path("results", "markers_MAST"))
markers_MAST2 = FindMarkers(integrated, ident.1=2, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST2, file.path("results", "markers_MAST", "cluster2-markers_MAST.csv"))
knitr::kable(head(markers_MAST2, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
pre-rRNA:CR45847 0 -1.6052927 0.435 0.771 0
CG13315 0 1.3570550 0.888 0.628 0
pre-rRNA:CR45845 0 -1.3662956 0.546 0.838 0
CG43349 0 1.3817990 0.643 0.296 0
pre-rRNA:CR45846 0 -1.8422695 0.354 0.625 0
MtnC 0 0.9374882 0.911 0.673 0
RpLP1 0 0.6907439 0.939 0.844 0
RpS26 0 0.6451404 0.942 0.857 0
Acbp5 0 1.3897104 0.870 0.672 0
RpS21 0 0.7083381 0.885 0.773 0
Cam 0 0.9412521 0.872 0.708 0
RpLP2 0 0.5753640 0.988 0.948 0
CG34330 0 0.8335108 0.868 0.654 0
28SrRNA:CR45844 0 -1.1858404 0.391 0.629 0
RpL41 0 0.6366086 0.990 0.972 0
RpL37A 0 0.6184733 0.965 0.914 0
RpL19 0 0.6111552 0.950 0.862 0
RpS30 0 0.6082303 0.935 0.857 0
RpS27 0 0.5947775 0.954 0.885 0
RpL13 0 0.6758232 0.871 0.746 0
RpL35 0 0.6034513 0.965 0.889 0
RpL23 0 0.5709788 0.965 0.920 0
RpS28b 0 0.6617437 0.903 0.800 0
18SrRNA:CR45838 0 -1.4476648 0.340 0.544 0
18SrRNA:CR41548 0 -1.4466951 0.340 0.544 0

Cluster 3

markers_MAST3 = FindMarkers(integrated, ident.1=3, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST3, file.path("results", "markers_MAST", "cluster3-markers_MAST.csv"))
knitr::kable(head(markers_MAST3, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
lncRNA:CR40469 0 2.302270 0.999 0.907 0
lncRNA:CR34335 0 2.055420 0.996 0.804 0
Tet 0 1.955220 0.607 0.167 0
sty 0 1.890713 0.616 0.261 0
cib 0 1.869054 0.716 0.303 0
E(spl)mbeta-HLH 0 1.832528 0.498 0.121 0
Df31 0 1.714096 0.662 0.234 0
E(spl)malpha-BFM 0 1.699416 0.438 0.115 0
hdc 0 1.682541 0.539 0.210 0
alphaTry 0 -3.590280 0.302 0.877 0
Jon65Aiv 0 -3.274725 0.089 0.748 0
esg 0 1.623873 0.627 0.165 0
betaTry 0 -3.749718 0.116 0.765 0
lncRNA:CR44898 0 1.181112 0.417 0.093 0
MtnC 0 -2.795806 0.108 0.738 0
mew 0 1.489866 0.524 0.196 0
robo2 0 1.338593 0.374 0.084 0
Vha16-1 0 -2.289525 0.243 0.837 0
CG30025 0 -3.694229 0.121 0.724 0
Jon99Ciii 0 -3.770872 0.230 0.779 0
Bace 0 -3.387552 0.118 0.728 0
CG13315 0 -2.205541 0.068 0.695 0
N 0 1.341854 0.367 0.102 0
Jon99Cii 0 -3.174153 0.155 0.699 0
Jon65Aiii 0 -2.990330 0.205 0.740 0

Cluster 4

markers_MAST4 = FindMarkers(integrated, ident.1=4, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST4, file.path("results", "markers_MAST", "cluster4-markers_MAST.csv"))
knitr::kable(head(markers_MAST4, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
Mal-A7 0 2.5066108 0.846 0.268 0
CG7542 0 1.9216820 0.917 0.277 0
Bace 0 1.6662209 1.000 0.657 0
Diedel3 0 1.4456951 0.986 0.462 0
CG5107 0 1.3861378 0.976 0.488 0
epsilonTry 0 1.3547922 0.967 0.434 0
deltaTry 0 1.3768449 0.974 0.504 0
CG30031 0 1.3478686 0.972 0.507 0
gammaTry 0 1.3478686 0.972 0.507 0
CG30025 0 1.3850072 1.000 0.652 0
betaTry 0 1.3489190 1.000 0.693 0
CG8834 0 1.2595427 0.873 0.261 0
Mal-A1 0 1.6094906 0.919 0.459 0
Cyp4e1 0 1.5399858 0.788 0.314 0
alphaTry 0 1.2351138 1.000 0.820 0
Pebp1 0 1.2475968 0.974 0.567 0
CG4377 0 1.3900148 0.820 0.302 0
Jon65Aiv 0 1.0712409 0.992 0.674 0
Mal-A8 0 1.4995523 0.663 0.198 0
CG4363 0 1.3258361 0.683 0.231 0
Jon65Aiii 0 0.8930227 0.991 0.676 0
tobi 0 1.5496700 0.744 0.289 0
MtnA 0 -2.7459049 0.410 0.731 0
LysD 0 1.1315302 0.605 0.207 0
CG12374 0 0.8567660 0.946 0.665 0

Cluster 5

markers_MAST5 = FindMarkers(integrated, ident.1=5, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST5, file.path("results", "markers_MAST", "cluster5-markers_MAST.csv"))
knitr::kable(head(markers_MAST5, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
LManVI 0 3.2826160 0.979 0.243 0
LManV 0 2.9089739 0.915 0.111 0
ninaD 0 2.6390118 0.871 0.126 0
CG15534 0 2.2950064 0.852 0.135 0
LManIII 0 2.2069193 0.697 0.052 0
LManII 0 2.1743168 0.869 0.141 0
CG31343 0 2.0399043 0.972 0.294 0
Mal-A4 0 1.9131132 0.813 0.161 0
iotaTry 0 1.8477332 0.836 0.170 0
LManI 0 1.7075755 0.760 0.109 0
CG7025 0 1.7006117 0.792 0.132 0
Jon99Cii 0 -2.6760089 0.608 0.660 0
CG14629 0 1.4508569 0.741 0.135 0
yip7 0 -2.6590889 0.577 0.544 0
Jon65Aiv 0 -2.9567791 0.605 0.704 0
CG15533 0 1.4710121 0.688 0.063 0
Pebp1 0 -2.5958549 0.586 0.598 0
Cyp6d5 0 1.6516791 0.891 0.332 0
Jon65Aiii 0 -2.9732678 0.612 0.705 0
lambdaTry 0 1.5033115 0.862 0.217 0
Npc2d 0 1.7057162 0.684 0.113 0
Rrp46 0 0.2793078 0.541 0.033 0
CG11911 0 1.2775764 0.940 0.346 0
asRNA:CR45281 0 1.1363389 0.746 0.158 0
Jon99Ciii 0 -3.2440617 0.649 0.743 0

Cluster 6

markers_MAST6 = FindMarkers(integrated, ident.1=6, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST6, file.path("results", "markers_MAST", "cluster6-markers_MAST.csv"))
knitr::kable(head(markers_MAST6, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
CG43774 0 3.1529264 0.807 0.207 0
MtnD 0 2.6718753 0.775 0.307 0
Arc1 0 2.4032135 0.911 0.374 0
MtnC 0 1.7872526 0.973 0.675 0
CG15423 0 1.9412040 0.737 0.166 0
mbl 0 2.1404858 0.850 0.388 0
thetaTry 0 1.9437544 0.623 0.136 0
CG34330 0 1.4639846 0.945 0.655 0
CG5773 0 1.8438525 0.352 0.028 0
CG5399 0 1.6595736 0.857 0.351 0
CG34301 0 2.0971596 0.504 0.156 0
Gagr 0 0.5444069 0.415 0.024 0
CG15422 0 1.7420240 0.730 0.245 0
CG8177 0 2.2074678 0.574 0.251 0
CG15127 0 1.7561668 0.354 0.076 0
Stat92E 0 1.7067355 0.631 0.212 0
MtnB 0 1.6650008 0.708 0.303 0
CG5770 0 1.6160867 0.315 0.023 0
MtnA 0 1.4680532 0.959 0.693 0
pre-rRNA:CR45847 0 -1.3870317 0.571 0.755 0
MtnE 0 1.6779141 0.642 0.249 0
28SrRNA:CR45844 0 -1.2991613 0.504 0.617 0
Tsp42Ec 0 1.5186541 0.562 0.187 0
CG9672 0 1.7261596 0.306 0.112 0
Adat1 0 1.0285194 0.640 0.249 0

Cluster 7

markers_MAST7 = FindMarkers(integrated, ident.1=7, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST7, file.path("results", "markers_MAST", "cluster7-markers_MAST.csv"))
knitr::kable(head(markers_MAST7, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
Jon65Aii 0 3.293820 0.814 0.210 0
Jon99Fi 0 3.266696 0.862 0.228 0
CG7916 0 3.257496 0.750 0.252 0
Jon25Biii 0 3.256755 0.960 0.411 0
Jon99Fii 0 3.236508 0.940 0.255 0
Jon65Ai 0 2.885554 0.858 0.142 0
CG7953 0 2.720265 0.800 0.203 0
CG17571 0 2.704814 0.936 0.257 0
CG10472 0 2.616301 1.000 0.405 0
mag 0 2.295322 0.929 0.172 0
CG3106 0 1.981872 0.638 0.041 0
CG8997 0 3.178884 0.801 0.263 0
Ag5r2 0 1.913701 0.658 0.116 0
CG11911 0 2.168190 0.933 0.348 0
CG17633 0 2.406228 0.811 0.222 0
PGRP-SC1b 0 0.844151 0.454 0.007 0
CG16749 0 1.553809 0.894 0.276 0
Acbp5 0 1.379220 0.998 0.670 0
CG8093 0 1.531794 0.670 0.092 0
Jon66Ci 0 2.028641 0.552 0.041 0
Jon44E 0 2.093519 0.517 0.065 0
CG3868 0 1.400003 0.902 0.330 0
Bace 0 -3.270406 0.466 0.694 0
CG31198 0 1.287501 0.885 0.290 0
asRNA:CR45281 0 1.095471 0.709 0.161 0

Cluster 8

markers_MAST8 = FindMarkers(integrated, ident.1=8, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST8, file.path("results", "markers_MAST", "cluster8-markers_MAST.csv"))
knitr::kable(head(markers_MAST8, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
CG5767 0 4.1823644 0.949 0.170 0
MtnB 0 3.4244227 0.947 0.297 0
CG30479 0 3.3674916 0.898 0.135 0
CG6277 0 3.3655225 0.646 0.112 0
CG8661 0 3.3205974 0.757 0.145 0
Vha100-4 0 2.9384083 0.664 0.047 0
CG31446 0 2.6010173 0.653 0.064 0
CG15423 0 2.5263314 0.895 0.165 0
Tsp42Ec 0 2.5020736 0.833 0.179 0
Adat1 0 2.4189735 0.889 0.242 0
Vha16-1 0 2.3347834 1.000 0.783 0
CG30480 0 2.2971498 0.802 0.056 0
Vha13 0 1.8668964 0.996 0.667 0
CAH1 0 2.2784895 0.777 0.132 0
CG17109 0 2.0909091 0.644 0.130 0
Vha55 0 2.3064912 0.958 0.416 0
CG5399 0 2.1605505 0.913 0.354 0
CG11192 0 0.3031397 0.474 0.005 0
CG31087 0 2.0721848 0.751 0.208 0
Vha14-1 0 1.9189858 0.900 0.299 0
Vha44 0 2.2445595 0.964 0.474 0
Vha68-2 0 2.1241612 0.931 0.341 0
CG15422 0 1.9328725 0.884 0.243 0
VhaAC39-1 0 2.2447769 0.913 0.328 0
CG6271 0 1.1021336 0.428 0.005 0

Cluster 9

markers_MAST9 = FindMarkers(integrated, ident.1=9, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST9, file.path("results", "markers_MAST", "cluster9-markers_MAST.csv"))
knitr::kable(head(markers_MAST9, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
Amy-p 0 3.7546156 0.926 0.413 0
CG6839 0 3.6581995 0.622 0.156 0
CG3819 0 3.3652770 0.657 0.126 0
CG11400 0 3.3465943 0.688 0.146 0
CG4928 0 3.3092892 0.678 0.197 0
CG14125 0 4.9471279 0.563 0.197 0
CG15818 0 3.2587424 0.492 0.108 0
Amy-d 0 3.4545206 0.718 0.233 0
CG15043 0 3.0710478 0.642 0.159 0
Mco1 0 2.4300353 0.424 0.061 0
CG13324 0 2.1388592 0.355 0.133 0
CG7589 0 2.2166688 0.386 0.075 0
Ag5r 0 2.4361222 0.307 0.066 0
CG13323 0 2.4128346 0.612 0.282 0
Adhr 0 1.6133663 0.541 0.202 0
CG42486 0 1.0222549 0.310 0.015 0
CG14933 0 1.4466203 0.368 0.045 0
CG18223 0 0.6894215 0.231 0.006 0
Adh 0 1.5147031 0.503 0.189 0
CG10405 0 1.0750603 0.246 0.027 0
lncRNA:CR40469 0 0.6702035 0.995 0.911 0
CG14439 0 0.7704120 0.241 0.010 0
CG42729 0 0.2725671 0.165 0.001 0
isoQC 0 1.4823857 0.358 0.228 0
CG9701 0 0.5324651 0.208 0.007 0

Cluster 10

markers_MAST10 = FindMarkers(integrated, ident.1=10, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST10, file.path("results", "markers_MAST", "cluster10-markers_MAST.csv"))
knitr::kable(head(markers_MAST10, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
18SrRNA-Psi:CR45861 0 1.9929792 0.774 0.425 0
18SrRNA:CR41548 0 1.8092135 0.917 0.516 0
18SrRNA:CR45838 0 1.8088839 0.917 0.516 0
18SrRNA:CR45841 0 1.8323186 0.902 0.506 0
pre-rRNA:CR45846 0 1.8543423 0.967 0.592 0
18SrRNA-Psi:CR41602 0 1.9522780 0.714 0.415 0
pre-rRNA:CR45856 0 1.9342579 0.842 0.497 0
pre-rRNA:CR45847 0 1.6122745 0.994 0.737 0
pre-rRNA:CR45845 0 1.5562434 0.979 0.810 0
28SrRNA-Psi:CR45851 0 1.3256885 0.688 0.531 0
28SrRNA-Psi:CR41609 0 1.3964064 0.560 0.459 0
28SrRNA:CR45844 0 1.1993126 0.744 0.607 0
28SrRNA-Psi:CR40741 0 1.3577011 0.565 0.448 0
28SrRNA:CR45837 0 0.7395038 0.411 0.486 0
CG12374 0 0.3401385 0.943 0.678 0
28SrRNA-Psi:CR45855 0 0.4141037 0.098 0.243 0
14-3-3epsilon 0 -0.6633088 0.158 0.488 0
Jon99Cii 0 0.5646288 0.893 0.650 0
Hr4 0 0.5679335 0.321 0.432 0
rg 0 0.3604389 0.137 0.328 0
Bace 0 0.2581049 0.923 0.674 0
His2Av 0 -0.3896197 0.062 0.330 0
Jon99Ciii 0 0.2742590 0.938 0.731 0
CG30025 0 0.2913271 0.899 0.671 0
CG42322 0 -0.2620372 0.182 0.471 0

Cluster 11

markers_MAST11 = FindMarkers(integrated, ident.1=11, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST11, file.path("results", "markers_MAST", "cluster11-markers_MAST.csv"))
knitr::kable(head(markers_MAST11, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
Acbp3 0 2.7989831 0.978 0.298 0
CG32473 0 2.6193899 0.739 0.143 0
Gs2 0 2.5571275 0.978 0.369 0
CG15254 0 2.5365965 0.796 0.138 0
zetaTry 0 2.2494411 0.873 0.180 0
CG31266 0 2.3784007 0.672 0.081 0
NAAT1 0 1.9975076 0.669 0.058 0
CG8774 0 2.0324823 0.768 0.148 0
CG13492 0 1.8982504 0.876 0.217 0
CG31343 0 1.8074765 0.984 0.310 0
CG4053 0 1.8113458 0.573 0.030 0
Pdk 0 1.7664866 0.640 0.155 0
CG31267 0 2.2780699 0.564 0.030 0
CG31265 0 1.8975320 0.580 0.052 0
Gdh 0 1.7188263 0.768 0.192 0
Pepck2 0 1.7380777 0.653 0.096 0
CG9673 0 1.7049865 0.895 0.234 0
CG31269 0 1.3435526 0.433 0.007 0
CG31091 0 1.9780214 0.503 0.031 0
CG17475 0 1.9069711 0.411 0.005 0
CG32379 0 1.2844718 0.427 0.008 0
CG5958 0 1.7785181 0.726 0.177 0
CBP 0 0.2538569 0.344 0.001 0
Agpat4 0 1.5228590 0.761 0.194 0
CG10116 0 1.7645785 0.812 0.219 0

Cluster 12

markers_MAST12 = FindMarkers(integrated, ident.1=12, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST12, file.path("results", "markers_MAST", "cluster12-markers_MAST.csv"))
knitr::kable(head(markers_MAST12, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
CCHa1 0 2.954678 0.603 0.119 0
IA-2 0 2.647187 0.993 0.403 0
nrv3 0 2.247964 0.909 0.242 0
Npc2b 0 2.621671 0.593 0.121 0
CG42541 0 2.090760 0.547 0.137 0
7B2 0 2.096422 0.915 0.231 0
w 0 2.291362 0.759 0.206 0
CCHa2 0 4.295649 0.407 0.249 0
Galphao 0 1.910813 0.691 0.175 0
pros 0 2.015971 0.769 0.181 0
Pal2 0 2.079188 0.782 0.147 0
CG42613 0 1.576964 0.450 0.043 0
CG30183 0 1.964119 0.736 0.130 0
Phm 0 1.793405 0.866 0.234 0
heph 0 1.849458 0.664 0.104 0
CG17646 0 1.608020 0.453 0.094 0
CG7191 0 2.818244 0.332 0.074 0
Tk 0 3.468051 0.557 0.218 0
Mrp4 0 1.856973 0.528 0.050 0
AstC 0 3.204312 0.524 0.302 0
amon 0 1.943279 0.638 0.127 0
Nrx-1 0 1.644528 0.570 0.095 0
CG32547 0 1.679896 0.463 0.056 0
unc-13-4A 0 2.278254 0.632 0.172 0
Hsp22 0 1.719189 0.720 0.342 0

Cluster 13

markers_MAST13 = FindMarkers(integrated, ident.1=13, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST13, file.path("results", "markers_MAST", "cluster13-markers_MAST.csv"))
knitr::kable(head(markers_MAST13, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
CG42825 0 2.2269327 0.761 0.116 0
CG14499 0 3.0172332 0.764 0.075 0
lectin-37Da 0 2.7082242 0.757 0.081 0
Mur29B 0 2.5662759 0.804 0.184 0
CG33926 0 2.1060157 0.950 0.228 0
CG10912 0 2.2182013 0.907 0.203 0
CG32368 0 2.1692952 0.854 0.168 0
asRNA:CR44192 0 1.3426450 0.635 0.150 0
Npc2e 0 2.1977119 0.691 0.154 0
CG10911 0 1.5862509 0.963 0.401 0
CG31086 0 1.5386260 0.804 0.311 0
CG8353 0 1.6389178 0.535 0.155 0
CG31323 0 1.7430702 0.837 0.335 0
CG11912 0 2.4572814 0.412 0.086 0
CG10943 0 1.3910626 0.498 0.151 0
CG18493 0 1.4654841 0.754 0.253 0
CG9568 0 1.6351462 0.688 0.224 0
CG15255 0 1.3555371 0.771 0.267 0
lectin-37Db 0 1.3335771 0.688 0.226 0
RpL3 0 -0.9937845 0.718 0.896 0
CG13482 0 1.4957212 0.625 0.217 0
Cht4 0 0.7473134 0.309 0.016 0
Nazo 0 1.1559679 0.518 0.128 0
CG8950 0 0.3780558 0.455 0.065 0
CG14500 0 0.8257770 0.518 0.090 0

Cluster 14

markers_MAST14 = FindMarkers(integrated, ident.1=14, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST14, file.path("results", "markers_MAST", "cluster14-markers_MAST.csv"))
knitr::kable(head(markers_MAST14, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
NPF 0 5.1479016 0.983 0.290 0
chrb 0 3.5003289 0.841 0.281 0
IA-2 0 3.2135054 1.000 0.403 0
7B2 0 3.0448420 0.924 0.231 0
Phm 0 2.6654247 0.890 0.234 0
Tk 0 2.4481634 0.857 0.210 0
svr 0 2.5674166 0.867 0.243 0
Ih 0 2.4344804 0.661 0.097 0
CG15312 0 1.9534549 0.591 0.071 0
CG30183 0 2.1819966 0.764 0.129 0
nrv3 0 2.2510264 0.827 0.245 0
cpo 0 2.4463486 0.761 0.275 0
Hsp23 0 2.8419301 0.781 0.318 0
esg 0 2.2221630 0.761 0.184 0
Pal2 0 2.1567229 0.748 0.149 0
Hsp22 0 1.8522887 0.777 0.341 0
Sh 0 1.6991290 0.492 0.042 0
Hk 0 1.8939359 0.495 0.026 0
Ldh 0 2.2306744 0.728 0.177 0
CG46385 0 2.2108067 0.831 0.325 0
unc-13-4A 0 1.7157612 0.757 0.169 0
AstA-R2 0 1.3280768 0.425 0.014 0
Hsp70Aa 0 1.7015774 0.724 0.352 0
Or49b 0 0.8860271 0.342 0.002 0
Hsp70Ab 0 1.6730323 0.714 0.350 0

Cluster 15

markers_MAST15 = FindMarkers(integrated, ident.1=15, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST15, file.path("results", "markers_MAST", "cluster15-markers_MAST.csv"))
knitr::kable(head(markers_MAST15, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
Orcokinin 0 5.3271710 1.000 0.292 0
AstC 0 3.3331658 0.996 0.293 0
CG4587 0 2.9631960 0.743 0.008 0
CG14989 0 3.7994074 0.871 0.190 0
CG13135 0 0.5614222 0.680 0.002 0
svr 0 3.0093252 0.979 0.244 0
CG33639 0 0.2836585 0.676 0.002 0
tap 0 1.9801462 0.851 0.030 0
RYa-R 0 0.3823188 0.593 0.000 0
unc-13-4A 0 2.9465008 0.963 0.167 0
Pvf3 0 0.3481559 0.676 0.005 0
CG30340 0 0.8985034 0.614 0.002 0
blot 0 0.6021470 0.680 0.009 0
Rgk3 0 2.1701331 0.768 0.022 0
Gad1 0 0.2640746 0.573 0.001 0
CG30116 0 1.2765806 0.714 0.013 0
dac 0 1.2441220 0.710 0.013 0
CG34458 0 0.7437293 0.685 0.010 0
Slob 0 0.7153288 0.685 0.011 0
brp 0 1.3592246 0.739 0.021 0
CG15203 0 0.8878283 0.544 0.001 0
CG15894 0 0.3864234 0.676 0.011 0
Frq1 0 0.5909662 0.614 0.009 0
CG3955 0 0.3744050 0.660 0.011 0
CG4341 0 0.5783358 0.685 0.014 0

Cluster 16

markers_MAST16 = FindMarkers(integrated, ident.1=16, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST16, file.path("results", "markers_MAST", "cluster16-markers_MAST.csv"))
knitr::kable(head(markers_MAST16, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
EbpIII 0 4.834704 0.920 0.223 0
spidey 0 3.329530 0.734 0.216 0
CG31522 0 3.284502 0.688 0.226 0
to 0 3.587573 0.726 0.064 0
Adhr 0 3.095023 0.819 0.200 0
Adh 0 3.083127 0.776 0.187 0
CG1124 0 3.623984 0.620 0.110 0
CG30197 0 3.532219 0.684 0.094 0
emp 0 2.438345 0.523 0.089 0
wat 0 3.231272 0.532 0.046 0
CG10237 0 2.613266 0.536 0.048 0
ple 0 2.881654 0.553 0.030 0
Msr-110 0 2.666640 0.612 0.116 0
CG31523 0 2.969051 0.658 0.148 0
ADPS 0 2.858588 0.473 0.050 0
Acbp2 0 2.549511 0.903 0.364 0
vir-1 0 2.554232 0.494 0.033 0
CG3097 0 2.313336 0.485 0.065 0
ATPCL 0 2.223677 0.527 0.170 0
Men 0 2.158755 0.515 0.191 0
CG8306 0 2.509803 0.494 0.060 0
CG43134 0 4.113872 0.354 0.116 0
Ace 0 1.772902 0.430 0.135 0
Cpr67B 0 2.731262 0.435 0.042 0
CG4660 0 2.353659 0.376 0.042 0

Cluster 17

markers_MAST17 = FindMarkers(integrated, ident.1=17, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST17, file.path("results", "markers_MAST", "cluster17-markers_MAST.csv"))
knitr::kable(head(markers_MAST17, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
Skp2 0 5.3000419 1.000 0.389 0
CG34324 0 5.2816359 0.987 0.277 0
CG14645 0 5.1663525 1.000 0.282 0
CG34220 0 5.0936717 0.987 0.420 0
CG3906 0 5.0871549 0.991 0.217 0
Muc68D 0 5.0305282 0.933 0.207 0
CG11672 0 3.8428054 0.871 0.037 0
CG4783 0 3.2040502 0.956 0.160 0
Pgant4 0 2.4506775 0.662 0.014 0
Idgf4 0 2.4551618 0.782 0.072 0
CG30026 0 1.7591050 0.511 0.003 0
CG43673 0 2.9267167 0.467 0.020 0
opm 0 1.6150495 0.462 0.147 0
CG9988 0 1.2403834 0.351 0.001 0
CG31077 0 2.0244394 0.298 0.002 0
sgl 0 1.4856293 0.387 0.187 0
Cys 0 1.5106883 0.507 0.172 0
Pgant8 0 0.8396921 0.271 0.002 0
CG13810 0 0.9988078 0.329 0.008 0
CG45061 0 1.3525408 0.342 0.026 0
Gfat1 0 1.3953112 0.360 0.097 0
NijC 0 1.1420618 0.271 0.018 0
Ctl2 0 1.3710576 0.387 0.104 0
Hexo1 0 1.4742258 0.489 0.142 0
CG30025 0 -3.1632394 0.116 0.690 0

Cluster 18

markers_MAST18 = FindMarkers(integrated, ident.1=18, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST18, file.path("results", "markers_MAST", "cluster18-markers_MAST.csv"))
knitr::kable(head(markers_MAST18, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
Npc2f 0 3.6101816 0.988 0.150 0
thetaTry 0 3.0983381 0.988 0.149 0
CG4830 0 2.6994819 0.789 0.017 0
CG18404 0 4.0870292 0.845 0.161 0
CG4563 0 2.4116003 0.770 0.061 0
Peritrophin-15a 0 2.6028515 0.950 0.184 0
CG9682 0 2.5655588 0.727 0.027 0
yip7 0 2.3389780 1.000 0.539 0
Jon65Aiv 0 2.1546865 1.000 0.694 0
Jon65Aiii 0 2.0314714 1.000 0.695 0
CG4734 0 1.8411515 0.826 0.132 0
CG30031 0 1.3828164 1.000 0.535 0
gammaTry 0 1.3828164 1.000 0.535 0
Bace 0 1.8774625 0.994 0.678 0
deltaTry 0 1.3582514 1.000 0.532 0
CG34040 0 0.5422801 0.640 0.082 0
LysB 0 1.3202083 0.826 0.197 0
28SrRNA-Psi:CR45853 0 0.3463674 0.422 0.025 0
Pebp1 0 1.5259964 0.944 0.592 0
CG10477 0 1.4562131 0.665 0.131 0
Pgcl 0 0.4049199 0.553 0.069 0
lncRNA:CR34335 0 -2.0779373 0.652 0.821 0
LysD 0 1.1953521 0.776 0.228 0
CG30025 0 1.0353096 1.000 0.673 0
LManII 0 0.9659594 0.708 0.172 0

Cluster 19

markers_MAST19 = FindMarkers(integrated, ident.1=19, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST19, file.path("results", "markers_MAST", "cluster19-markers_MAST.csv"))
knitr::kable(head(markers_MAST19, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
CG6409 0 4.803699 0.966 0.145 0
CG42656 0 4.436801 0.910 0.082 0
IM33 0 3.965089 0.831 0.081 0
CG1143 0 2.846481 0.584 0.052 0
Irk2 0 2.881903 0.584 0.057 0
Ggamma30A 0 2.425736 0.584 0.061 0
tws 0 2.662732 0.640 0.170 0
CG14949 0 1.958402 0.315 0.027 0
CG16704 0 2.138345 0.360 0.030 0
Nha2 0 2.023512 0.382 0.042 0
Gdh 0 2.581828 0.685 0.205 0
CG9993 0 1.838703 0.371 0.012 0
Vha13 0 1.570255 0.978 0.679 0
CG18473 0 1.749039 0.337 0.018 0
Msr-110 0 2.274386 0.596 0.123 0
Gdap1 0 1.572230 0.270 0.149 0
Nop17l 0 1.449145 0.326 0.154 0
GLS 0 1.632297 0.180 0.046 0
lncRNA:CR40469 0 1.381869 1.000 0.913 0
bru1 0 1.433266 0.247 0.014 0
CG14933 0 1.864400 0.461 0.054 0
CG17999 0 1.509553 0.303 0.011 0
sesB 0 1.222895 0.978 0.775 0
Vha16-1 0 1.169516 0.978 0.791 0
CG7365 0 1.059547 0.225 0.015 0

Cluster 20

markers_MAST20 = FindMarkers(integrated, ident.1=20, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST20, file.path("results", "markers_MAST", "cluster20-markers_MAST.csv"))
knitr::kable(head(markers_MAST20, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
CG13285 0 3.8504561 0.963 0.020 0
CG44013 0 2.9005776 0.778 0.020 0
Msr-110 0 3.1788920 0.963 0.123 0
Idgf6 0 2.6238014 0.778 0.043 0
CG43394 0 2.6414936 0.685 0.028 0
GstD11 0 2.0589475 0.537 0.001 0
CG10096 0 2.8104611 0.667 0.012 0
CG5065 0 2.2587751 0.741 0.095 0
Idgf4 0 2.3408393 0.815 0.083 0
ELOVL 0 1.5613344 0.370 0.087 0
bond 0 2.0227798 0.481 0.016 0
Skp2 0 2.3129974 0.593 0.401 0
CG9336 0 1.7430299 0.500 0.041 0
CG17839 0 1.6639478 0.426 0.030 0
CG34220 0 2.2185041 0.556 0.432 0
yellow-d 0 1.4443748 0.333 0.002 0
Cys 0 1.8452262 0.611 0.177 0
trn 0 1.5692802 0.407 0.090 0
Gasp 0 1.5855129 0.389 0.008 0
Muc68D 0 1.8044905 0.519 0.221 0
pyr 0 1.1379466 0.296 0.205 0
CG34324 0 1.9788417 0.500 0.291 0
Baldspot 0 1.9471257 0.537 0.252 0
Snp 0 0.9650331 0.167 0.131 0
CG14645 0 2.2891081 0.593 0.296 0

Cluster 21

markers_MAST21 = FindMarkers(integrated, ident.1=21, test.use="MAST") %>%
  as.data.frame() %>%
  tibble::rownames_to_column("gene") %>%
  arrange(p_val_adj)
readr::write_csv(markers_MAST21, file.path("results", "markers_MAST", "cluster21-markers_MAST.csv"))
knitr::kable(head(markers_MAST21, 25))
gene p_val avg_logFC pct.1 pct.2 p_val_adj
Mur18B 0 3.763457 0.583 0.080 0
Alp2 0 3.229087 0.556 0.029 0
CG6726 0 3.195073 0.611 0.109 0
Alp4 0 4.009171 0.639 0.055 0
CG3168 0 3.775266 0.528 0.073 0
Spat 0 1.457579 0.278 0.026 0
Oatp58Dc 0 1.726650 0.222 0.015 0
CG2680 0 1.974845 0.361 0.026 0
CG10513 0 2.711481 0.361 0.037 0
CG7882 0 2.203185 0.167 0.016 0
CG10226 0 1.227024 0.250 0.032 0
Smvt 0 2.016667 0.306 0.015 0
CG10553 0 1.840484 0.333 0.016 0
CG7084 0 2.335300 0.361 0.091 0
CG6733 0 1.778015 0.306 0.028 0
CG10560 0 1.147329 0.139 0.020 0
Irk3 0 1.842469 0.333 0.018 0
CG31373 0 1.888380 0.333 0.112 0
CG10514 0 2.743189 0.333 0.009 0
Tps1 0 1.900322 0.222 0.009 0
HDAC6 0 2.516269 0.500 0.161 0
UGP 0 2.228555 0.500 0.160 0
CG13309 0 2.391852 0.333 0.008 0
CG14292 0 3.042500 0.417 0.098 0
CG11426 0 1.755369 0.306 0.109 0
saveRDS(integrated, file.path("results","integrated.rds"))